Under consideration for publication in J. Fluid Mech. 



1 



Dynamics of reorientations and reversals of 
large scale flow in Rayleigh-Benard 
convection 

P. K. MISHRA 1 , A. K. DE 2 , 
M. K. VERMA'and V. ESWARAN 3 

department of Physics, Indian Institute of Technology, Kanpur 208 016, India 

2 Department of Mechanical Engineering, Indian Institute of Technology Guwahati 781039, 

India 

3 Department of Mechanical Engineering, Indian Institute of Technology, Kanpur 208016, India 

(Received 27 August 2010) 

We present a numerical study of the reversals and reorientations of the large scale cir- 
culation (LSC) of convective fluid in a cylindrical container of aspect ratio one. We take 
Prandtl number to be 0.7 and Rayleigh numbers in the range from 6 x 10 5 to 3 x 10 7 . 
It is observed that the reversals of the LSC are induced by its reorientation along the 
azimuthal direction, which arc quantified using the phases of the first Fourier mode of the 
vertical velocity measured near the lateral surface in the mid plane. During a "complete 
reversal", the above phase changes by around 180° leading to reversals of the vertical 
velocity at all the probes. On the contrary, the vertical velocity reverses only at some 
of the probes during a "partial reversal" with phase change other than 180°. Numeri- 
cally we observe rotation-led and cessation-led reorientations, in agreement with earlier 
experimental results. The ratio of the amplitude of the second Fourier mode and the 
first Fourier mode rises sharply during the cessation-led reorientations. This observation 
is consistent with the quadrupolar dominant temperature profile observed during the 
cessations. We also observe reorientations involving double cessation. 

Key Words: Rayleigh-Benard Convection, Convective Turbulence, Direct Numerical 
Simulation. 



1. Introduction 

Turbulent convection is ubiquitous in nature and in many engineering applications. 
Rayleigh-Benard Convection (RBC) in which fluid confined between two plates is heated 
from below and cooled on the top is an idealized yet an important paradigm to understand 
convective turbulence. The dynamics of RBC is governed by the two non-dimensional 
parameters: the Rayleigh number R = aATd 3 g/iyn and the Prandtl number P = v/k, 
where d is the vertical height of the container, g is the acceleration due to gravity, 
AT is the temperature difference between the bottom and top plates, and a, k, and v 
are the thermal heat expansion coefficient, thermal diffusivity, and kinematic viscosity, 
re spectively, of the fluid. 



Krishnamurti fc Howard! (|198ll ) performed experiments on water (P ~ 7.0) and silicon 
oil (P ~ 860) and observed coherent roll structur es, also known as "lar ge scale circula- 



tion" (LSC), in the turbulent regime. Subsequently. ICastaing et all ( 19891 ) ascertained the 
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existence of LSC in Helium (P ~ 0.65-1.5) contained in a cylindrical container. They pro- 
posed that coherent large scale structures exist statistically only above a certain Rayleigh 
Number (i? ~ 10 s ). T hey also observed a low frequency peak in the power spectrum of 
the temperature field. Xi. Lam, fc Xia ( 2004 ) studied the onset of large-scale coherent 
mean flow in RBC using shadowgraph and particle image velocimetry techniques and 
s howed that LSC is a result of the org anization of plume motion. 



Cioni. Ciliberto. fc Sommerial (|1997f ) performed RBC experiments on mercury (0.021 < 



P < 0.026) and placed several thermistors along the azimuth of the cylinder. They de- 
duced the presence of global circulation from the dipolar temperature distribution mea- 
sured by the probes. The temperature fluctuations (after subtracting the mean) switched 
sign randomly in their experiment. Since the warmer fluid ascends from one side, and the 
cooler fluid descends from the other side of the apparatus. ICioni. Ciliberto. fc Sommeria 
( 19971 ) deduced that the vertical velocity would also exhibit random "reversals" in phase 

with the temperature fluctuations. This feature of convection has been s tudied extensively 

using theoretical, exp erimental, and computational tools (see reviews bv lAhlers. Grossmann. fc Lohse 



2009t lKadanoflll200lh . In this paper we perform computational investigation of the re- 



v ersal dynamics of LSC in a cylindrica l geometry. 



Cioni. Ciliberto. fc Sommerial (|1997l ) computed the first Fourier mode of the measured 



temperature field. They observed that the amplitude of the first Fourier mode never 
vanishes, but the phase of the Four ier mode is highly variable. A phase change of ir cor- 
responds to the reversal of the flow. Cioni. Ciliberto. fc Sommeria ( 19971 ) also observed a 
low-frequency peak in the energy spectrum that corresp onds t o the circulation frequency 
of the large-scale flow inside the cylinder. Tsuji et al\ (|2005l ) performed an experiment 
with mercury contained in cylindrical container for aspect ra tios 1/2, 1, and 2 and found 
that the low frequency peak is abse n t for aspect ratio 1/2. iNiemela et al. ( 2001 ). and 
Sreenivasan. Bershadskii. fc Niemela (120021) reported reve rsals in their convect i on ex- 



periments on Helium. Brown. Nikolaenko. fc Ahlersl ( 2005 ) and Brown fc Ahlerj ( 2006f ) 
performed similar experiments on water (P = 4.38) and measured the temperature field 
in the bulk using probes placed at l/4th, 1/2 and 3/4th height of the cylinder. They ob- 
served that the plane containing LSC exhibits diffusive and sometimes quick significant 
drift in the azimuthal direction with the angular change (|A0|) ranging from small values 
to 7r, i.e., 7r/4 < |A0| < 7r. This is called the "reorientation" of LSC. The reorientation of 
the flow can occur in two ways: (a) through rotation of the circulation plane without any 
major reduction of the circulation strength; (b) through " cessation" of the circula tion, 
followed by a restart in a randomly chosen new direction. iBrown fc Ahler 1 (|2006l ) also 
computed the amplitude of the first Fourier mode, and found it to be non zero for rotation- 
led re orientations and close to zero during cessation-led reorientations. IXi. Zhou, fc Xia 
(2006) measured azimuthal motion of the LSC using particle image velocimetry and stud- 
ied reversals of LS C. They ob s erved "double cessation" in th eir experiments. In a similar 
set of experiments. IXi fe Xial (f2007h . IXi fc Xial (|2008ah . and lXi fe Xial (|2008tJ ) observed 
st rong dependence of the azimuthal reorientations on the aspect ratio of the cylinder. 
Qiu fc Tond ( 2001 ) performed RBC experiments on water and computed temperature 



correlation functions. They reported a transition from ra ndom chaotic states to a corre- 
lated turbulent states at around R c ~ 5 x 10 7 . However. ISano. Wu. fc Libchaberl (jl989h 
report strongly correlated large scale flow beyond R c ~ 4 x 10 7 for helium gas (P = 0.7). 

Va rious models have been proposed to under stand the dynami cal behaviour of the 
LSC. ISreenivasan. Bershadskii. fc Niemela (|2002Tl and|Benzj| (l2005h proposed a stochas- 
tic model in which the reversal of LSC was explained as nois e induced switching between 
two meta-stable states. lAraujo. Grossmann. fc Lohsel (|2005h attempted to explain the ir- 
regular cessation and subsequent reversal of the LSC using the force and thermal balance 
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on a s ingle plume mod e lled b y co upled nonlinear equation s related to the Lorenz equa- 
tions. Brown fc Ahlerd (|2007t ) and iBrown fc Ahlersl (|2008l ) proposed a stochastic model 
with the strength and the azimuthal orientation of the LSC being determined by two 
stochastic ordinary differential equations; this model appear s to e xplain the experimen- 
tal observations of cessation and reorientation. Villermaus| (|l995 ) proposed a model to 
explain the low frequency oscillations in the temperature signal in the bulk. His model 
assumes that the modes of the boundary layers interact through the slow moving large 
scale circulation. 

There are only a small number of c omputational studies on the reorientation or reversal 
of LSC. IStringano fc Verziccol (|2006t) simulated convection in air (P = 0.7) in a cylinder 
with an aspect ratio of 1/2 and observed a sing le roll breaking into two counter rotating 
rolls stacks vertically. Benzi & VerziccJ (|2008 ) performed a simulation of the fluid at 
R = 6 x 10 5 with white noise adde d to the heat equ a tion a nd studied the statistical 
behaviour of the observed reversals. Breuer fc Hansen (2009) studied RBC for infinite 
Prandtl number in a 2D box and observed reversals of LSC for very high R (R = 10 9 ). 
There are others 2D numerical e xperiments in a box that report that t he reversals of 
the LSC is due to cessations (see Hansen. Yuen, fc Kroenin j 1990l 19921) or due to the 
chaotic movement of rolls perpendicular to the roll axis (see lPaul et a/j|2010h . 

In this paper we study the dynamics of reversals of convective structures using direct 
numerical simulation (DNS) of turbulent RBC for P = 0.7 in a cylinder of aspect ra- 
tio one. We apply conducting boundary conditions on the top and bottom plates and 
adiabatic boundary conditions on the lateral walls of the container. No-slip boundary 
condition is applied on all the surfaces of the container. Our simulations indicate a LSC 
in the container. Furthermore, we observe rotation-led and cessation-led reorientations 
of LSC similar to those observed in experiments. We compute the amplitudes and phases 
of the Fourier modes of the vertical velocity measured near the lateral wall in the mid 
plane, and establish that the dynamics of the LSC can be captured quite well by the low 
wavenumber Fourier modes. 

The outline of the paper is as follows. We present our numerical method in Section 2. 
The results pertaining to observation of LSC and their reorientations are presented in 
Section 3. We conclude our results in Section 4. 



2. Numerical Method 

We numerically simulate the convective flow in a cylinder. The Boussinesq approxima- 
tion is assumed for the buoyancy in the fluid. The relevant non-dimesionalized dynamical 
equations for the fluid are 

5 + (u-V)u = -W + Tz+\/^V 2 u, (2.1) 



dt y ' V R 

dT 
~dt 



^- + (u • V)T = ^=V 2 T, (2.2) 



Vu = (2.3) 

where u = (u x , u y ,u z ) is the velocity field, T is the temperature field, a is the deviation of 
pressure from the conduction state, R = ag(AT)d 3 /isk is the Rayleigh number, P = v / 'k 
is the Prandtl number, and z is the buoyancy direction. Here v and k are the kinematic 
viscosity and thermal diffusivity respectively, d is the vertical height of the container, 
and AT is the temperature difference between the bottom and top plates. For the non- 
dimensionalization we have used d as the length scale, ^/a(AT)gd (free fall velocity) as 
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Table 1. N r , Ng, and N z are the number of grids along the radial, azimuthal, and vertical direc- 
tions of the cylindrical container; A m i n and A max are the minimum and maximum mean-grid 
sizes; lmax = max(2nr/Ng); rjh is the Kolmogorov length scale calculated using the expres- 
sion rjh — tt(P 2 /R * Nu) 1 / 4 ; Nbl is the number of points inside the thermal boundary layer; 
JVu(comp) is the Nusselt number obtained from the simulation; (e") (= ^(|Vu| 2 )) and (ejT) 
(= k(|VT| 2 )) are the numerically calculated viscous and thermal dissipation rates respectively; 
(ea) (= v a (Nu - l)RP~ 2 /d 4 ) and (e^) (= K,(AT) 2 Nu/d 2 ) are the analytical values of viscous 
and thermal dissipation rates respectively. The last two columns show (e")/(e") and (e^)/(t^). 



the velocity scale, and AT as the temperature scale. Consequently, d/^/agATd is the 
time scale of our simulation. The aspect ratio of the container is taken to be one. We 
confine our study to P = 0.7 which is a typical Prandtl number for air. 

The above non-dimensionalized equations (|2.HI2.3[) are solved numerically for a cylin- 
drical geometry using a finite-difference scheme. The convective parts of the equations are 
discrctized in cylindrical coordinates using Tar n and Webb's fourth order central explicit 
scheme with enhanced spectral resolution fsee lTam &; Webblll993l ). The diffusive part is 
discretized using the second-order central-difference scheme. For the time advancement, 
we use the second-order Adam-Bashforth scheme for the nonlinear terms, and the Crank- 
Nicholson scheme for the diffusive terms. We perform simulations for Rayleigh numbers 
R = 6 x 10 5 , 8 x 10 6 , 2 x 10 7 , and 3 x 10 7 . Comparison with experimental results show 
that these Rayleigh numbers are near the threshold of strong turbulence regimes (see 
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The cylinder volume is discretized into variable grids with finer resolution near the 
boundary layers. Since the boundary layers significantly affect the dynamics of LSC in 
convective turbulence, it is nec essary to resolve the regions nea r the top and bottom 
plates, and the lateral walls (see lStevens. Verzicco. fc LohsJ2010h . In our simulations we 
choose uniform grids along the azimuthal direction, and non-uniform on es along the radia l 
and vertical directions. For the grid spacing, the Grotzbach condition (|Grotzbachlll983h 
is used, according to which the mean grid size A = (rAOArAz) 1 ' 3 should be smaller 
than the Kolmogorov and thermal diffusion length scales. Note that the Kolmogorov 
length scale (rjh) is estimated using the formula rjh — tt(P 2 /(R * Nu)) 1 / 4 , where Nu is 
the Nusselt number. We also calculate the width of the thermal boundary layer using the 
formula 5t ~ 2NH anc ^ ensure that the number of grid points inside the boundary layer 
(Nbl) should be greater than 3 to 5. Number of grid points inside the thermal boundary 
layer is given in Table [TJ and they satisfy the above condition. 

Our simulations satisfy the Grotzbach condition as the minimum of the mean grid size 
(A rn i„) is smaller than r/h (see Table[T]). Note however that the maximum grid size l m ax = 
max(2TTr /Ng) is sometimes larger than rjh- To validate our code, we compare the numeri- 
cally computed Nusselt numbers at two different grids (one more resolved than the other) 
and find that the Nusselt numbers for the two different grids are quite close. We also cal- 
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Figure 1. (a) Some of the velocity and temperature probes placed inside the cylinder. For most 
of our discussions in the paper we use the probes at z = 0.5, r = 0.48 and 9 = jV/4 with j — : 7 
shown in the figure, (b) Temperature profile in a vertical plane 9 = 37r/4 for R = 2 x 10 7 on 
a (100 x 120 x 201) grid. A hot plume ascends from the right wall and a cold plume descends 
from the left wall confirming the presence of large scale structure, (c) Temperature profile in a 
horizontal section z = 0.5 for R = 2 x 10 7 on a (100 x 120 x 201) grid. 



culate the kinetic energy dissipation rate (e") (= ^(|Vu| 2 )) and the thermal dissipation 
rate (e^T) (= k(|VT| 2 )) using the numerical data and compare them with their theoretical 
estimates (e£) (=v 3 (Nu ~ l)RP- 2 /d 4 ) and (e£) (=K(AT) 2 Nu/d 2 ) (IShraiman. fc Siggial 
1990t) . The ratios (e£) /{el) and (ej) /(e^) are listed in the Tabled] Clearly the numerical 



kinetic dissipation rate is quite close to its analytical counterpart, with the maximum 
difference at 12%. For the thermal dissipation, the numerical value always seems to be 
higher than the analytical value, with the maximum difference around 35%. However 
note that the theoretical estimates of (e") and (e;T) are functions of the Nusselt number 
for which we substitute its numerically computed value. Ou r numerical results on dissi- 
pation rates appear to be less accurate compared to those bv lStevens. Verzicco. fc Lohse 
( 2010( ). Also, our numerical thermal dissipation rates are always l arger than the theoreti- 
cal estimate in contrast to those by Stevens. Verzicco. fc Lohse ( 2010l ) who report lower 
values compared to the corresponding theoretical estimate. The above discrepancies are 
probably due to the diff erence in the time integration sche mes and the spatial derivative 
schemes used by us and Stevens. Verzicco. fc Lohse] ( 201dt ). 

We use fixed time-step dt which is listed in Tabic [T] Due to the computational com- 
plexity, the longest computer run for R = 2 x 10 7 took approximately 45 days on 32 
cores of CHAOS cluster of IIT Kanpur. To ensure that the initial conditions do not af- 
fect our final results, we start every run with a conductive state modulated with random 
noise. We validate our numerical code by comparing our results with the past numerical 
results. For example, our numerical results of Nusselt number Nu for different Rayleigh 
numbers fit with a relatio n Nu ~ 0.143i? 0,297 whic h is i n good agreement with th e 
earlier numerical results of Stringano fc Verzicco! (I2006TI and I Verzicco fc Camussil (|l999f ). 



and the experimental observations oflNiemela et al. ( 2000). Similar to some of the ear 



Brown. Nikolaenko. fc Ahlersl 



her experiments (see ICionh Ciliberto. fc Sommerialll997t ^_ 
l2005t iBrown fc Ahlersl 1200a IXi fc Xiall2008bl ). we place "probes" in the bulk and in the 
boundary layers of the cylinder to record the time series of the velocity and temperature 
fields. The probes are located at 9 = 6j = jir/A (j = : 7), the vertical heights of 
z = 0.02, 0.25, 0.5, 0.75, 0.98, and the radial distances of r = 0.15, 0.3, 0.48 from the axis. 
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Table 2. N r , Ne, and N z are the number of grids along the radial, azimuthal, and vertical 
directions of the cylindrical container; u^ ms is the rms value of the vertical speed of the flow; 
and T c ddy = 2d/u z lns is the circulation time of large scale flow. 



The bottom and top probes are inside the thermal boundary layers. A schematic diagram 
of the probes in the middle plane is depicted in Fig.[TJa). The nondimcnsional rms values 
of the velocity m™ s , and the eddy turn-over time (T e ddy = 2d/it™ s ) for various i?'s are 
listed in Table[2] u™ s has been computed by taking the time average of ("z) 2 ) 1 ^ 
for the eight probes placed at z = 0.5 and r = 0.48. The eddy turnover time ranges from 
26 to 10 free fall times for the range of Rayleigh numbers studied in the present paper; its 
va lue decreases with the increase of Rayleigh number consistent with the earlier results 
of iQiu fc Toriel (l200ll) . 

In the next section we will study the properties of the temperature and velocity time 
series, and relate them to the earlier experimental results on reversals of LSC. 



3. Numerical results 

3.1. Large scale circulation 

The time series measured by the probes of Fig. [IJa) carries signature of the LSC that 
will be described below. Later in the section we will relate the reversals of the vertical 
velocity field with the reorientation of the LSC. In Figs. [TJb) and QJc) we display the 
temperature profiles for R = 2 x 10 7 for a vertical plane [6 = 3tt/A) and the middle 
horizontal plane (z = 0.5) respectively. As a supplementary material, we also provide 
two short movies depicting the flow behaviour for the above two sections. These figures 
and movies clearly indicate the presence of convective structures in the flow. 

In Figs.[^a,b,c,d) we plot the normalized autocorrelation g c (T) = (ST(t)6T(t + T)) /erf* 
at 6 = 0,r = 0.48, z = 0.5_for R = 6 x 10 5 ,8 x 10 6 ,2 x 10 7 , and 3 x 10 7 resp_ectively. 
Here 5T(t) = T(t) - f (T is the mean temperature) and erf = T{t) 2 - (T) 2 , and 
time is measured in the units of eddy turnover time. The autocorrelation functions in- 
dicate two time-scales in the system. Variations at the shorter time-scale, shown in the 
insets, are due to the statistical return of the convective flow after one eddy turnover 
time or less (see Table [TJ. For R = 6 x 10 5 and 8 x 10 6 , the oscillations are somewhat 
irregular. However for R = 2 x 10 7 and 3 x 10 7 , the oscillations in the insets are quite 
regular, and the oscillation time period of the auto-correlation function is around one 
eddy turnover time. This is also evident from the plot of the cross-correlation function 
of temperature between two azimuthally opposite probes placed at 9 — and 6 = tt, 
r = 0.48, and z = 0.5 (see Fig. [3]). The cross-correlation function also shows oscillations 
with approximately one eddy turno v er tim e as the time period, in general a g reement with 
the ob serva tions of Castaing et all ( 1989f ). Sano. Wu. fc Libchaber ( 1989I ). Qiu fc Tong 
(j200lh . and lXi. Zhou, fc Xial (|2006T ). The above correlations are related to the peak in the 
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Figure 2. Plot of the autocorrelation function g c (r) = {5T(t)5T(t + T))/a^ vs. delay time r (in 
units of large eddy turnover time) for a probe at r — 0.48, 9 = 0, and z = 0.5: (a)i? = 6 x 10 5 
on a (33 x 49 x 97) grid, (b)R = 8 x 10 6 on a (75 x 96 x 145) grid, (c)R = 2 x 10 7 on a 
(100 x 120 x 201) grid, and (d)R = 3 x 10 7 on a (100 x 120 x 201) grid. The insets show 
oscillations on eddy turnover time scale. 




Figure 3. For R = 2 x 10 7 on a (100 x 120 x 201) grid, the cross-correlation function between 
the temperature measured at two azimuthally opposite probes at 9 = and 9 = n in the mid 
plane at r = 0.48. The inset shows oscillations on eddy turnover time scale. 



power spectrum at frequency corresponding to 1/T e ddy The peak in the spectrum how- 
ever tends to be overshadowed by noise. Consequently the correlation function appears 
to be a good tool to analyze these oscillations. Note that the autocorrelation studies 
tend to beco me strongly periodic when LSC direc tion is locked, e.g., by a small tilt of 
the cylinder ( Brown. Funfschilling. fc Ahlersll2007 ). 

In addition to the above, we observe variations in the above auto-correlation and cross- 
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t (free fall time) 

Figure 4. For R = 6 x 10 5 on a (33 x 49 x 97) grid, the time series of the vertical velocity 
measured by the probes at z = 0.5 (mid plane), r = 0.48: (a) 9 = and 9 = ir, (b) 9 — 7r/4 
and 9 = 57r/4, (c) 9 = tt/2 and 9 = 3tt/2, and (d) = 3tt/4 and 9 = 7vr/4. The time series of 
the phase of the first Fourier mode of the vertical velocity (8i) is shown in subfigure (e). Time 
is measured in the units of the free fall time dj \/agATd. The arrows (inside boxed region) in 
(e) indicate a partial reorientation (Si ~ 135 degree) near t ~ 2700, and a complete reversal 
(<5i = 180) near t ~ 3900. 

correlation functions at much larger time-scales. For R = 6 x 10 5 (see Fig. [2fa)), there is 
a significant decrease in correlations till r ~ 2ZT e( idyi after which correlations rise again 
to reach a maximum value at t ~ h%T e ^ y . These variations possibly correspond to the 
reorientations of the LSC. Similar features are observed for R = 8 x 10 6 (see Fig. H^b)). 
For R = 2 x 10 7 and 3 x 10 7 , the auto-correlation functions decrease with time (see 
Figs. [2lc,d)). Incidentally, the cross-correlation function for R = 2 x 10 7 shows minima 
at t = and then at t w 20T ee [d y , and a maximum at t ~ 10T ec idy, which possibly 
indicates reversals of the flow at r w l0T e ddy, and reoccurrence at t ps 20T e ddy- We 
could not perform longer statistics since the computer simulations of RBC flows for large 
Rayleigh numbers are computationally very expensive. The auto-correlation and cross- 
correlation function studies are inconclusive due to averaging, yet, the above observations 
indicate with some certainty that the LSC exist in RBC. The reversal time could not 
be deduced from the present correlation functions. The presence of LSC becomes more 
apparent when we will study the time series of the velocity and temperature fields, a 
topic of next subsection. 

3.2. Measure of reorientations of LSC 

The time series of the vertical velocity field for R = 6 x 10 5 , 8 x 10 6 , 2 x 10 7 recorded by 
the probes at dj = jn/4 (j = 0: 7), r = 0.48, and z = 0.5 are shown in Figs. |H [SJ and [6] 
respectively. Here time is measured in the units of the free fall time [dj -Jag(AT)d). The 
figures clearly show that vertical velocity at probes 8j and (dj +ir) mod 27r are clearly anti- 
correlated. These observations indicate the presence of a large scale convective structure 
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Figure 5. For R = 8 x 10 6 on a (75 x 96 x 145) grid, the time series of the vertical velocity 
measured by the probes at z = 0.5 (mid plane) and r = 0.48, and the phase of the first Fourier 
mode of the vertical velocity. Details of the figures are same as Fig. 0] A complete reversal of 
the flow occurs in the boxed region. 



(see Figs. [Hb,c)). As shown in Figs. |4]|6]the local mean value (in time) of the vertical 
velocity changes sign, and this feature is called "reversal" of LSC in literature. 

In the following discussions we will show that the above mentioned reversals are con- 
nected to the "reorientations" of LSC. To quantify the reorientations, we Fourier trans- 
form the vertical velocity field measured at r = 0.48 and z — 0.5. Note that these probes 
are near the lateral wall in the mid plane. The velocity signal at a given probe can be 
expressed in terms of its Fourier transform as 



j(t) = ^2 Ukexp(ik9j) 

fc=-4 
4 

= 2J \u k \ cxp (ikd j + 5 k ) 

k—i 

^2|w fe | cos(k6j + 6 k ), 



fe=-4 



(3.1) 
(3.2) 
(3-3) 



fc=i 



where 6j = jir/ifj = 0:7), and 8 k is the phase of the fc-th Fourier mode. Note that the 



reality condition u_ 



fit is used to derive Eq. 



In our simulations we observe that 



the first Fourier mod e (fc = 1) dominates over the other mo des e xcept possibly dur i ng the 
reversals. Note that Cioni. Ciliberto. fc Sommeria (1997), and Brown fc AhlersI ( 20061 ) 
have analyzed the properties of u\(t) and Si(t) only using their experimental data. In 
later part of our discussion we will compute the amplitudes of the Fourier modes using 

the vertical velocity fields m easured by probes at z = 0.5 and r = 0.48.^ 

For quantitative analysis, Brown. Nikolaenko. fc AhlersI ( 2005 ) and Brown fc AhlersI 
(200G) proposed that the phase of the first Fourier mode <5i can be used as an approx- 
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Figure 6. For R = 2 x 10 7 on a (100 x 120 x 201) grid, the time series of the vertical velocity 
measured by the probes at z = 0.5 (mid plane) and r = 0.48, and the phase of the first Fourier 
mode of the vertical velocity. Details of the figures are same as Fig. 0] A complete reversal of 
the flow occurs in the boxed region. 



imate measure for the orientation of the plane of LSC. Using this criteria, they defined 
"reorientation of the LSC" as a sudden and significant rotation of the above plane. They 
used two selection criteria for the reorientation of LSC: (a) the magnitude of angular 
change in orientation (|A<$i|) should be greater 7r/4; (b) the azimuthal rotation rate 
{\S\ |) should be greater than ir / (5T e dd y ) , where T e dd v is the eddy turn over time. We 
follow the same criteria for the selection of reorientation events in our simulation. 

In Figs. |4]|6] we plot the phase of the first Fourier mode (6±), which is a measure of 
the plane of LSC. Note that the discontinuities from 180 degree to —180 degree in the Si 
diagrams are not reorientations; they simply indicate jitters near 180 degree. In all the 
time series of the vertical velocity and the phase of the first Fourier mode, we observe that 
the mean value of the vertical velocity changes sign but not necessarily simultaneously 
for all the probes. For example in Fig. [4j the vertical velocity u z measured by all the 
probes reverse sign near t ~ 3900 (the right boxed region of Fig. |4}. This reorientation 
is "complete reversal" of LSC, and it corresponds to the change in 6± by around 7T, i.e., 
AS i « 7T. Near t ~ 2700 (the left boxed region of Fig. HJ however u z changes sign for all 
the probes except at 9 = and 7r. This kind of reorientation will be termed as a "partial 
reversal", and it corresponds to AS\ ^ tt. Near t ~ 2700, A<5i « 135°. The time series 
shown in Figs. |4][6] exhibit several partial and complete reversals. For example, we observe 
complete reversals for R = 8 x 10 6 and 2 x 10 7 in the boxed regions shown in the figures. 
A careful analysis of the Fourier modes reveals that the complete and partial reversals are 
intimately related to the change in orientations of the convective structures by A^i . We 
take the real space u z data before and after the reversals for the three cases marked with 
arrows in the 5\ time series of Figs. |4] and [H In Figs. [7Ja,c,e) we illustrate the azimuthal 
profile of vertical velocity at z = 0.5 and r = 0.48 before and after the reorientations 
as blue and red curves respectively. We calculate A8\, the change in the phase of the 
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Figure 7. Azimuthal profile of the vertical velocity at z = 0.5 and r = 0.48 before (blue curve) 
and after the reorientation (red curve): (a) R = 2 x 10 7 at t — 695 and 750; (c)R = 6 x 10 5 
at t = 2548 and 2850; (e) R = 6 X 10 5 at t = 3400 and 4050. During the reorientations the 
phase of the first Fourier mode changes by ASi = 180, 135 and 180 degrees respectively. For the 
reoriented time series we subtract ASi from the positive k modes, and add ASi to the negative 
k modes, and construct a modified velocity profile (green curve). As seen in the subfigures (b), 
(d), and (f), the reconstructed velocity profiles match quite well with the profiles before the 
reorientation. 



first Fourier mode, during the reorientations of the structures. If the convective structure 
rotates by an angle A<5i during a reorientation, then, according to Eq. (|3.3p , we can cancel 
the effects of this reorientation by subtracting kAdi from the phases of all the positive 
k modes, and by adding fcA^i to all the negative k modes of the data recorded after 
the reorientation (Note that the k = mode is left unaltered). We perform the above 
exercise on the Fourier modes of the data sets recorded after the specified reversals or 
reorientations. Subsequently we compute the velocity field from the modified Fourier 
coefficients. The modified velocity field, shown as green curves in Figs.[3Jb,d,f), matches 
quite well with the corresponding velocity profiles before the reorientations. Thus we 
show that the convective structures essentially rotate by A5i during the reorientations. 
The ab ove features of our simulations are in general agreement with the experim e ntal 



result s by Cioni. Ciliberto. fc Sommerial ( 1997 ) , Brown. Nikolaenko. k, Ahlers ( 2005 ) , Brown fc Ahlerd 
( 20061 ) . and lXi. Zhou, fc Xial ( 2006| ). Note that the numbers of reorientations in the con- 
vection experiments are much larger than those observed in numerical simulations. In 
the next subsections we analyze the detailed dynamics of the reorientation of the LSC 
using numerical data. 



3.3. Dynamics of reorientations 

Careful analysis of the phases and amplitudes of the Fourier modes provide impor- 
tant clues about the dynamics of the reorientations of LSC and the reversals of the 
vertical velocity. The amplitudes of the Fourier modes vary significantly over time. 
We observe that during some of the reversals or reorientations, the amplitude of the 
first Fourier mode almost vanishes; these reorientations are termed as "cessation led". 
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Figure 8. Time series of the phase of the first Fourier mode (8i) (in units of 2ir), amplitude 
of the first Fourier mode (|«i|), and amplitude of the second Fourier mode (\u2\ ) during the 
rotation-led reversals, (a), (b), (c), and (d) show the time series of 8i, \ui\, |«2|, and |u2|/|tii| 
respectively for R = 6 x 10 s on a (33 x 49 x 97) grid in which the reversal takes place near 
t = 3950. (e), (f), (g), and (h) show the time series of t5i, \v,2\, and |u2|/|mi| respectively 
for R = 2 x 10 7 on a (100 x 120 x 201) grid with a reorientation event near t = 700. 



In the other reorientations that are termed as "rotation led" the Fourier modes con- 
tin ue to fluctuate around their mean v alues. These kinds of r eversa ls were reported 
bv iBrown. Nikolaenko. fc Ahlersl (|2005l ). iBrown fc Ahlersl (|2006l ). and IXi. Zhou, fc Xia 
( 20061 ) based on their RBC experiment. In the present paper we probe the dynamics of 
reversals using numerical simulations. In Fig.[S]wc plot the phases and amplitudes of the 
first two Fourier modes during two rotation-led reorientations. For R = 6x 10 5 Fig.[5Ja) 
exhibits the time series of the phase (5i/2ir, i.e., in units of 2ir). Figs. [8jb,c) show the 
amplitudes of the first and second Fourier modes (|mi|, |w2|) respectively, while Fig.^d) 
shows their ratio |u2|/|t2i|. As shown in the figures, Si changes by approximately 7r near 
t ps 3900, but and \u2\ continue to fluctuate around their average values. However, 
always dominates over \u2\, as evident from the plot of |u2|/|mi| shown in Fig. [5Jd). 
Figures [8]Je,f,g,h) exhibit similar features for R = 2 x 10 7 . Note however that tends 
to have a small dip during the reorientation, but A|mi| is much smaller compared to the 
cessation-led reversals to be described later. Note that both the above events lead to 
"complete reversals" since w it. 

In the other set of reversals or reorientations, we observe that the amplitude of the first 
Fourier mode (|mi|) t ends to vanish during the even t . The s e set of reversals were te rmed 
as "cessation led " bv | Brown. Nikolaenko. fc Ahlersl (120051 ). [Brown fc Ahlersl (|2006l ). and 



Xi. Zhou, fc Xial (120061 ). We illustrate the cessation-led reversals using the following ex- 



amples. 

Figures a) and [DJb) exhibit the time series plots of r5i and |ui| respectively for 
R = 8 x 10 6 . During the cessation at t <~ 1175, A5i ~ it and tends to zero. The 
amplitude of the second Fourier mode |u2| however increases slightly above its average 
value during this event (see Fig.[9jc)). The fluctuations in |{ti| and |u2| however tend to 
hide the above features. Therefore we use |m2|/|mi| to amplify the decrease in and 
increase in \v,2\ to be able to identify the cessation- led reversals clearly. In |9jd) we plot 
|w2|/|ui| that exhibits a sharp peak during the event. This is an example of cessation-led 
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Figure 9. Time series of Si, \ui\, \u2\, and |u2|/|ui| indicating cessation- led reorientation near 
t ~ 1175 for R = 8 x 10 6 (75 x 96 x 145). The strength of |«i| decreases significantly, while 
|w2|/|wi| increases during this event. 
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Figure 10. Time series of Si, \ui\, \u2\, and J«2|/|ui| indicating cessation-led partial reorienta- 
tions near t ~ 6625 and 6710 for R = 6 x 10^ (33 x 49 x 97 grid). For both the events Ar5i arc 
approximately 72°. 



complete reversal since A<5i ~ n. Figurc fTUI shows the corresponding quantities for another 
set of cessation- led reorientations for R = 6 x 10 5 . As evident from the time series, near 
t - 6620, AS 1 - (3.2 - 2) * 2tt w 72°. Later at t - 6720, Aft - (2 - 1.2) * 2tt w -72°. 
These two reorientations of Fig. [10] arc partial ones, and the LSC comes back to its 
original configuration after the two reorientations (A<5i = (3.2 — 1.2) * 2tt). During both 
these events |u2|/|ui| has a sharp peak. 

In our si mulations we also obser ve reorientations involving "double cessations" , first 
observed bv lXi. Zhou, fc Xial(l200ri in their experiments. During these events, the phase 
Si changes in two stages in quick succession, first by 9i and then second by 02 ■ Therefore 
the net AS± ~ 61+62- During double cessation, |tii| vanishes on two occasions separated 
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Figure 11. Time series plot of Si, \ui\, \u 2 \, and |«2|/|mi| indicating double cessation for 
R = 6 x 10 J (33 x 49 x 97 grid). First cessation occurs at t ~ 2083 and second one occurs 
at t ~ 2088. During both the events becomes weak and \u2\ dominates, (d) shows the two 
spikes in |u2|/|wi| during the events. 



by a small time gap (within an eddy turnover time). \u2\ tends to increase during these 
times, hence |zt2|/|ui| exhibits two peaks within a short time interval. Two independent 
double cessation events have been illustrated in Fig. HTT a) and Fig.[T27a) (for R = 6x 10 5 ). 
In Fig. HTT a) we observe the first cessation at t pa 2083 followed by the second cessation at 
t pa 2088, and 6\ pa — 62 ~ ir. In Fig. IT2T a) the first cessation occurs at t pa 2464 followed 
by the second cessation at t pa 2466 with 6\ pa — Q 2 pa 0.87T. Since the net change in the 
phase 5\ is approximately zero, the final configuration of LSC is similar to its original 
configuration. Note that one eddy turn time for the corresponding run is approximately 
26 free fall time (see Table 2). Hence the two double cessations occur within 0.2 and 
0. 08 eddy turnover time respectively, which is consistent with the experimental result 
of IXi. Zhou, fc Xial (|2006l ). 

For fi = 6x 10 5 we observe two events that are similar to the above mentioned double 
cessation events (see Fig. ITBT a.b.c.d)). Near t pa 4540, 8\ changes twice, 9\ pa 0.6-7T and 
82 pa — 0.67T with net A<$i pa 0. These two events however occur very close to each other, 
and we observe only one peak peak for |u2|/|ui|, rather than any double peaks indicative 
of double cessation. Another cessation-led reorientation occurs near t sa 4270 that resem- 
bles the above double-cessation event except one major difference. Here Q\ sa 0.87T and 
62 pa — 0.6-7T, hence LS C reorients by a net ang le after the event. Note that for the double 
cessation reported by XL Zhou, k, Xial ( 20061 ) Q\ sa —Q 2 and the net change in S\ is zero. 
In our numerical simulations we find such events. In addition we also observe a double ces- 
sation where 6\ ^ —02, and the LSC reorients by a finite angle after the event. Figure H4l 
illustrates the distribution of A<5i of LSC observed in our numerical simulations. We ob- 
serve three complete reorientations (ASi ~ n) among which two are rotation-led while one 
of them is cessation-led. The number of observed partial reorientations are rather large 
(total 19). Among the partial reorientation events, 15 of them are rotation-led, while 4 of 
them are cessation-led. We also observe 4 double cessations among whom three of them 
have AS 1 pa while the fourth one has ASi sa 0.27T. The above observations indicate that 
cessation- led eve nts are rarer compared to the rotatio n-led ones, an observation c onsis- 
terrt with those of lBrown. Nikolaenko. fc Ahlerd (|2005l ) and lBrown fc Ahlersl (|2006l ). The 
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Figure 12. Time series plot of Si, \iti\, \u2\, and |u2|/|ui| indicating double cessation for 
R = 6 x 10 5 (33 x 49 x 97 grid). First cessation at t ~ 2464 is followed by the second at 
t ~ 2466. During both the events |ui| becomes weak and \ii2\ dominates, (d) shows the two 
spikes in |w2|/|«i| during the events. 
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Figure 13. For R — 6 x 10 time series plot of 5i, |tta|, and |-u2|/|«i| depicting events 
similar to double cessations (33 x 49 x 97 grid). A<5i ~ 72° near t ~ 4270, while ASi w near 
t ~ 4540. During both the events |wi| — > but the ratio |w2|/|wi| increases significantly. 



number of reorientations observ ed in our simulations are far f ewer c o mpared to those 
observed in the experiments of Brown. Nikolaenko. fc Ahlers J 2005 \j_ Brown fc Ahlersl 



( 2006t ). ICioni. Ciliberto. fc Sommerial (|l997t ). and IXi. Zhou, fc Xial (l2006l). Consequently 
we are una ble to perform stat i stics o f reorientations similar to those by Brown. Nikolaenko. fc Ahlersl 
(|2005l ) and lBrown fc Ahlersl (|2006h . 

In our numerical simulations we observe that the higher Fourier modes play inter- 
esting role during reorientations. These kind of investigations were missing in earlier 
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Figure 14. The distribution of the change in the phase of the fi r st Fo urier mode (|Ar5i|) of 
LSCduring reorientations. Criteria of lBrown. Nikolaenko. fc Ahlersl (|2005f ) and lBrown fc Ah lcrs 
( 200(1) have been used for identifying reorientations of LSC. 



experimental works. We find that during the cessation-led reorientations |u2| tends to 
become significant when the |wi| tends to vanish. This feature is evident from the sig- 
nificant increase in the magni t ude o f fel /l&il during the cessat ion-led reorientations 
Brown. Nikolaenko. fc Ahlersl (|2005l ) and lBrown fc Ahlersl (|2006l ) interpret such events 



as cessation of the circulation, followed by a restart in a randomly chosen new direc- 
tion. Our simulation results are consistent with the above interpretation, yet another 
scenario is also possible. During the cessation, the circulation structure corresponding 
to iii (dipolar in a horizontal plane) becomes very weak, and the quadrupole structure 
corresponding to iii becomes dominant. After the cessation, the quadrupole structure 
disappears and dipolar structure reappears with a shift in the azimuthal direction. These 
features have been illustrated in Fig. [15] in which we show the horizontal profile of tem- 
perature in z — 0.5 plane before, during, and after two cessation-led reorientations for 
R = 6 x 10 5 . Figures [TST a.b.c) represent the temperature snapshots near a double ces- 
sation event (near t = 2100 of Fig. ITTj) . The system starts from ui dominant state with 
approximate profile as cos(# + Si) of Fig. ITST a). During the cessation, the system profile 
appears as Fig. 115( b) which has iii as the most prominent Fourier mode with profile as 
cos(26* + 62)- After the double cessation, the system return to iii dominant structure as 
evident from Fig. I15f c) . Similar features are observed for another cessation- led reorienta- 
tion near t ps 4500 for R = 6 x 10 5 shown in Figs. ITBT d.e.f). Note that in the rotation-led 
reorientations, the dipolar structure continues to be dominant during the reorientation 
itself. 

In this section we present some of the complex dynamics of the reorientations of LSC. 
The reversals of the vertical velocity have been shown to be intimately related to the 
reorientations of LSC . We observe rotation- led and cessa tion-l ed reorientations , whic h 
were first reported bv lBrown. Nikolaenko. fc Ahlers (l2005h and iBrown fc Ahlersl (120061) . 



Xi 



Xiall2007l 



2008b). 



and observed later in other experiments (|Xi. Zhou, fc Xiall2006 
In our simulations we also find double cessation, previously observed bv lXi. Zhou, fc Xia 
( 20061) . A new feature of our simulation is a double cessation in which the orientation 
of the LSC after the event is different from the original orientation. In our analysis we 
study higher Fourier modes, notably U2, which was missing in earlier work on reversals. 
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Figure 15. Temperature snapshots at the mid plane (z = 0.5) during cessation-led reorien- 
tations near t ~ 2087 and 4539 for R — 6 x 10 J . Dipolar structures are dominant before and 
after the reorientations, but quadrupolar structure is prominent during the reorientations (at 
t = 2087 and 4539). 



We propose a new interpretation of the cessation-led reorientations, and argue that the 
LSC transforms from dipolar (ui dominant) to quadrupolar (v,2 dominant), and then 
back to dipolar (ui dominant) structure during this event. 

4. Summary and Conclusions 

In summary, we numerically compute turbulent convective flows for a cylindrical geom- 
etry and study the characteristics of the LSC, in particular those related to flow reversals. 
The first Fourier mode along the azimuthal direction is the most dominant mode; the 
phase of this mode (<5i) is used as a measure of the orientation of LSC. Our numerical 
results and earlier experimental results indicate that the LSC occasionally reorients itself 
by rotating along the azimuthal direction by an arbitrary angle. When the reorienta- 
tion angle is around 7r, the vertical velocity at all probes change sign signaling complete 
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reversal of the flow. A partial reversal is observed when the reorientation angle is less 
than 7r. We therefore argue that the primary cause of the flow reversal is the azimuthal 
reorientation of the LSC. 

In our numerical simulations we obse rve two kinds of reorientations: (a) rotation-led, 

and (b ) cess ation-led, earlier observe d bv lBrown. Nikolaenko. fc Ahlersl ( 20051 ). Brown fc Ahlersl 
( 2006| ). and Xi. Zhou, fc Xia ( 20061 ). The rotation- led reversals involve rotation of the 

LSC without appreciable reduction in circulation strength, i.e., the ampli tude of |mi| re- 

mains finit e. In the cessation- l ed rev ersals, — > 0, which is interpreted bv lBrown. Nikolaenko. fc Ahlersl 
; 20051 ) and B rown fc Ahlersl (|2006l ) as cessation of circulation followed by a restart in a 
randomly chosen new direction. During this event, our numerical simulations reveal that 
the ratio of the amplitude of the second Fourier mode and the first Fourier mode, i.e., 
1^*2 1 / |^i | , increases significantly, and then it comes back to its original level. These prop- 
erties of the Fourier modes reveal that during the cessation, the LSC transforms from 
a dipolar-like structure to a quadrupolar-like structure, and then back to a dipolar-like 
structure. These features appear to have certain similarities with the reversals of magnetic 
field in dynamo (jWicht fc Olsonl 120041) . The role of higher Fourier modes for reversals 
and reorientations of LSC have been highlighted in our simulations for the first time. 

We also observe double cessation in our simulations. IXi. Zhou, fc Xia (200G) had ob- 
served in their experiments that the LSC returns to its original orientation (approxi- 
mately) after a double cessation. We however find in our numerical simulations that the 
change of the orientations of LSC during some double cessations are zero, while in some 
others, the changes are nonzero. 

Our numerical simulations reproduced many features ("rotation-led reorientations", 
"cessation- led reorientations", "double c e ssatio n" etc.) observed in convection exper - 
iments by Cioni. Ciliberto. fc Sommeria 



Brown fc Ahlersl (|2006f ) . IXi. Zhou, fc Xia 




Brown. Nikolaenko fc Ahlersl ([20051) . 



, and lNiemela et al\ (|200l[ ) for Rayleigh 



numbers greater than those in our simulations. Several researchers have argued that the 
reorientations, including reversals, of LSC occur only under strong turbulence regimes 
(say R > 10 s ). We however observe that the nature of convective flows and reorienta- 
tions are very similar for the range of Rayleigh numbers studied by us (R = 6 x 10 5 to 
R = 3 x 10 7 ). The lowest Rayleigh number R = 6 x 10 5 is probably in the weak turbu- 
lence regime where the LSC is not well organized, and R = 3 x 10 is in the lower end of 
strong turbulence regime. We could not carry out simulations for even larger Rayleigh 
numbers due to very expensive computational requirements for these runs. However we 
believe that reversals in high Rayleigh number regimes are highly likely to be similar to 
those presented in the simulations. 
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